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ABSTRACT 

Context. The Gaia catalogue will contain observations and physical parameters of a vast number of objects, including ultra-cool 
dwarf stars, which we define here as stars with a temperature below 2500 K. 

Aims. We aimed to assess the accuracy of the Gaia T^fj and log(g) estimates as derived with current models and observations. 
Methods. We assessed the validity of several inference techniques for deriving the physical parameters of ultra-cool dwarf stars: 
Gaussian processes, support vector machines, fc-nearest neighbours, kernel partial least squares and Bayesian estimation. In addition, 
we tested the potential benefits of data compression for improving robustness and speed. We used synthetic spectra derived from ultra- 
cool dwarf models to construct (train) the regression models. We derived the intrinsic uncertainties of the best inference models and 
assessed their validity by comparing the estimated parameters with the values derived in the bibliography for a sample of ultra-cool 
dwarf stars observed from the ground. 

Results. We estimated the total number of ultra-cool dwarfs per spectral subtype, and obtained values that can be summarised 
(in orders of magnitude) as 400000 objects in the M5-L0 range, 600 objects between LO and L5, 30 objects between L5 and 
TO, and 10 objects between TO and T8. A bright ultra-cool dwarf (with r^ff =2500 K and log(g)=3.5) will be detected by Gaia 
out to approximately 220 pc, while for Teff =1500 K (spectral type L5) and the same surface gravity, this maximum distance 
reduces to 10-20 pc. We found the cross-validation RMSE prediction error to be 10 K for regression models based on the 
A:-nearest neighbours and 62 K for Gaussian process models in the faintest limit (Gaia magnitude G=20). However these 
values correspond to the evaluation of the regression models with independent test sets of synthetic spectra of the same model 
families as used in the training phase (internal errors). For the ^-nearest neighbours model, this seems an overly optimistic 
error estimate due to the use of a dense grid of examples in the training set, together vnth a relatively high signal-to-noise ratio 
for the end-of -mission data. The RMSE of the prediction deduced from ground-based spectra of ultra-cool dwarfs simulated at 
the Gaia spectral range and resolution, and for a Gaia magnitude G=20 is 213 K and 266 K for the models based on fc-nearest 
neighbours and Gaussian process regression, respectively. These are total errors in the sense that they include the internal 
and external errors, with the latter caused by the inability of the synthetic spectral models (used for the construction of the 
regression models) to exactly reproduce the observed spectra, and by the large uncertainties in the current calibrations of 
spectral types and effective temperatures. We found maximum-likelihood methods (minimum fe-nearest neighbours, and 
Bayesian estimation vrith flat priors) to be biased in the LO-TO range in that they systematically assign a temperature around 
1700 K. Finally, the likelihood landscape is significantly multimodal in spectra with realistic noise. 

Key words. Methods: data analysis,Methods: statistical. Catalogues, brown dwarfs, Stars: fundamental parameters 



1. Introduction 

In this work we define an ultra-cool star as a star with an ef- 
fective temperature below 2500 K (spectral type M8). The goal 
of this paper is to assess the detectability of this type of ob- 
ject with the Gaia spacecraft. Gaia is a mission of the European 
Space Agency that will produce very accurate astrometry and 
parallaxes for a significant fraction of the galactic population (de 
Bruijne 2012), thus helping to considerably improve our knowl- 
edge of a plethora of astronomical topics, from stellar evolution 
to exoplanets. In particular, Gaia data will improve our under- 
standing of the nature of ultra-cool dwarfs by providing dis- 
tances and therefore luminosities for the nearest objects. The 
2MASS All-Sky survey (Cutri et al. 1996) began a new era in 



the discovery and characterisation of very-late spectral type stars 
and brown dwarfs, allowing the identification of two new types: 
L (Kirpatrick et al. 1999; Martin et al. 1997) and T (Kirpatrick 
et al. 1999; Burgasser et al. 2002), and paving the way for even 
cooler objects, the Y type (Burningham et al. 2008; Kirkpatrick 
et al. 2012). Subsequent surveys have discovered hundreds of 
cool objects, but the comprehensive understanding of their na- 
ture relies in modelling from internal structure to atmospheres. 
This can only be achieved with precise data, including accurate 
distances. Here, Gaia will truly play the role of a Rosetta stone. 

The work presented here was developed in the framework of 
the eighth coordination unit (CUB) of the Gaia Data Processing 
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and Analysis Consortium (DPAC ), which is devoted to deter- 
mining astrophysical parameters. Gaia is expected to detect and 
characterise one billion sources, and hence, automatic proce- 
dures for the reduction and processing of these data are essential. 

The DPAC consortium is, in charge of the design, develop- 
ment, and operation of this data processing and analysis chain 
aimed of producing the Gaia catalogue (intermediate releases 
and final catalogue) from the telemetry data (see Mignard et al. 
(2008) for a more detailed introduction to the DPAC). 

Since potential Gaia targets include very diff'erent astrophys- 
ical scenarios, from unresolved galaxies and quasars to aster- 
oids, specialised modules have been designed and implemented 
within the DPAC to characterise the various object types. 

The astrophysical parameters in CU8 are determined by var- 
ious modules integrated in the the Apsis pipeline. Apsis includes 
an initial classification into broad object categories (Bailer- Jones 
et al. 2008), modules for derivating stellar parameters from 
(amongst other observables) very low- and medium-resolution 
spectra (Liu et al. 2012; Recio-Blanco et al. 2006, respectively), 
and speciaUsed modules for characterising of unresolved galax- 
ies (Tsalmantza et al. 2012), quasars, or peculiar types of stars 
such as emission-line stars (Blomme et al. 2011), or cool stars. 

In particular, a specific module of the Gaia processing 
pipeline is devoted to characterising ultra-cool dwarfs (hereafter 
UCDs), which constitutes a regression problem in which the 
source parameters are estimated from observational data. In the 
next section, we describe this module and estimate the number 
of sources that will be detected as a function of spectral type. 
In Sect. 3 we briefly describe the statistical techniques explored 
in the search for an optimum model for the regression problem 
of determining the source parameters, together with the experi- 
ments carried out in order to select amongst them. In Sect. 4 we 
describe the results obtained by these models when they are ap- 
pUed to simulated Gaia spectra of well-known ultra-cool dwarfs 
observed with ground telescopes. These results provide a pre- 
launch approximation to the expected accuracy of the Gaia pa- 
rameter estimates (mainly for the efi'ective temperature). Finally, 
Sect. 5 summarises the main results of this work and describes 
the experiments that are currently being carried out to complete 
this study. 

2. Gaia sample of ultra-cool dwarfs 

2. 1 . Brief description of tiie Gaia capabilities 

The Gaia astrometric mission was approved by the European 
Space Agency in 2000 and the construction of the spacecraft and 
payload is on-going for a launch in mid 2013. Gaia will contin- 
uously scan the entire sky for five years, yielding positional and 
velocity measurements with the accuracies needed to produce 
a stereoscopic and kinematic census of about one biUion stars 
throughout our Galaxy and beyond. The stellar survey will be 
complete to Gaia magnitude G = 20 mag, with a precision of 24 
juas at magnitude V=15 for a solar-type star (G2V). Gaia will be 
equipped with two spectrophotometers operating in the 330-680 
nm range (blue photometer or BP) and in the 640-1000 imi range 
(red photometer or RP). Both spectrophotometers are based on 
a dispersive-prism approach, and the spectral resolutions are, in 
both cases, wavelength dependent. The RP has a varying resolu- 
tion from 7 nm pixel"' at 640 nm to 15 nm pixel"' at 1000 nm, 
while the BP photometer resolution reaches from 4 to 32 nm 



' http://www.rssd. esa. int/index.php?project= 
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pixel ' in its wavelength range. The details of the Gaia passband 
G and the photometric performances of the Gaia instruments are 
surmnarised in Jordi et al. (2010). 

2.2. Ultra-cool dwarfs with Gala 

The spectral energy distributions (SEDs) of ultra-cool dwarfs all 
peak in the infrared range, and we do not expect any significant 
flux in the BP range (see Fig. 1). Therefore, we will mainly be 
concerned with detecting and characterising of these stars using 
RP spectra. 

A rough estimate of the number of ultra-cool dwarfs that will 
be detected by Gaia per bin of spectral type can be obtained 
asuming a local volume density such as the one compiled by 
CabaUero et al. (2008). These authors compiled (or derived) lo- 
cal densities of late-type stars and brown dwarfs between spec- 
tral types M3 and T8 from the literature, to provide estimates 
of contamination rates by these objecs in deep photometric sur- 
veys that searched for substellar objects or high-redshift quasars. 
The overall shape of the local density as a function of spectral 
type shows three local maxima at M3, L5, and T8, and reaches 
a minimum of 0.22x10"^ objects per cubic parsec at spectral 
types TO-Tl. For each spectral type, we computed the distance at 
which a main-sequence object of the corresponding I-band abso- 
lute magnitude would reach the Gaia detection hmit of G=20 us- 
ing the Gaia object generator. This (hereafter GOG) is one of the 
three Gaia generators of simulated data that also include GASS 
(telemetry generator) and GIBIS (image generator). The GOG is 
a tool designed to obtain directly simulated catalogue and main 
database data for the Gaia satellite, passing through the entire 
mission data reduction chain (Robin et al. 2012). The outputs 
are astrometric, photometric, and spectroscopic epoch and final 
data. To simulate the main database data lifecycle, GOG uses 
error models whose formulas are coded using the current knowl- 
edge of the Gaia mission performances. In this work, we used 
GOG simulations of two synthetic libraries of ultra-cool dwarf 
spectra to derive the maximum distance at which an ultra-cool 
dwarf can be detected as a function of the I-band absolute mag- 
nitude, and other properties of the Gaia UCD sample. 

The first library is a composite of the AMES-Cond and 
AMES-Dusty models described in Allard et al. (2001). The va- 
lidity ranges for these models are Teff < 1400 K (AMES-Cond) 
and Jeff > 1700 K (AMES-Dusty). Therefore, there is a gap in 
the validity (not in the coverage) in efl'ective temperature be- 
tween 1400 and 1700 models. Models in the interregnum are 
available in both model famihes (and hence, no gap in cover- 
age exists), and wiU be used in this work to interpolate between 
the validity domains. The second Ubrary is the BT-Settl family 
of models (Allard et al. 2012), valid across the entire range of 
efTective temperatures. Figure 2 shows the evolutionary tracks 
for the BT-Settl library in the Tgs-log(g) space for a range of 
masses between 0.0005 and 1.4 Mq. Throughout this work we 
measure eff'ective temperatures in Kelvin and gravities in cm s"^. 
The simulation of the synthetic spectra is carried out in practice 
by the so-called coordination unit 2 (CU2) of the Gaia DPAC. 
In most of this work we concentrated the results obtained with 
the BT-Settl Ubrary of models which produced better fits to the 
observed spectra used for validation in Sect. 4. 

These two model libraries will also be used to define the 
mapping between the source parameters and Gaia observations 
described in section 3. Figure 1 shows a sample of spectra from 
the BT-Settl library together with the simulated Gaia RP spec- 
tra for a range of temperatures between 500 K and 2500 K. We 
obtained Gaia simulated data using the GOG. In the simulation 
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Fig. 1: Normalised sample spectra from the BT-Settl library. The top row shows simulated Gaia RP spectra of BT-Settl models for 
Teff = 500, 1000, 1500, 2000, and 2500 K. The vertical axis is proportional to the number of photons detected in each wavelength 
bin. The line colours reflect the various values of log(g) available in the library of models according to the colour scale on the right. 
The bottom row shows the original spectra with the same temperatures as in the top row, \og(g) = 5.0 and solar metallicity. 
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Fig. 2: Evolutionary tracks in the Tf.ff-log(g) space for ultra-cool 
dwarfs according to the BT-Settl library. Each line corresponds 
to a different mass in the range 0.0005-1 .4Mq as labelled below 
selected tracks. Filled circles represent individual models in the 
grid. These are coloured according to the decimal logarithm of 
the age in Gigayears as indicated by the colour scale on the right. 
Effective temperatures are measured in Kelvin and gravities in 
cm s"^. 



nitudes and effective temperatures of the BT-Settl model fam- 
ily, combined with the calibration of effective temperatures with 
spectral types by Stephens et al. (2009). We used the analytic 
formula based on optical spectral types of L dwarfs and infrared 
spectral types of the T dwarfs. This calibration (hereafter refered 
to as SLC calibration) is valid in the M6 to T8 range. This anal- 
ysis results in the values listed in Table 1 and illustrated in Fig. 
4(a). Table 1 lists the relationship between spectral type and ab- 
solute magnitude in the I-band given in Caballero et al. (2008) 
in column 2 for reference. This results in the estimated number 
of counts under the column header CountscBKQf,- The relation- 
ship between spectral type and absolute magnitude in the I-band 
implicit to the BT-Settl model family is included in column 3. 
This results in the estimated number of counts under the column 
header Count SBT-Setii- The SLC spectral type-effective temper- 
ature calibration (see Sect. 4) is included in column 6, and the 
(Gbp-Grp) colour index computed from the noiseless GOG sim- 
ulations of BT-Settl model atmospheres is included in column 7. 
Since the volume densities tabulated in Caballero et al. (2008) 
refer only to the main sequence, these estimates do not take into 
account the potential detection of low-gravity objects (i.e., es- 
sentially very young objects). The final increase in the expected 
number of counts is due to the steep increase in the volume den- 
sity for spectral types later than TO. In deriving these estimates 
we only used BT-Settl models of solar metallicity. The expected 
number of detections per apparent G magnitude is shown in Fig. 
4(b). 



of the libraries, CU2 takes into account the evolutionary tracks 
provided with the models. 

Figure 3 shows the maximum distances to main-sequence 
UCDs as a function of their I-band absolute magnitude. The 
expected number of detectable main-sequence objects in each 
spectral type bin can be estimated multiplying a volume density 
estimate (in our case, the one in Caballero et al. 2008, for spec- 
tral types between M3 and T8) by the volume of a sphere with 
a radius equal to the maximum distance at which an UCD cor- 
responding to the spectral type under consideration can be de- 
tected (assuming solar metallicity). In computing these expected 
number counts per spectral type bin, we need to define a rela- 
tionship between absolute I-band magnitude and spectral type. 
We used two such relations. The first one is included in Table 3 
of Caballero et al. (2008) (and reproduced in Table 1 for conve- 
nience). The second relationship is derived from the I-band mag- 



2.3. Selection criteria and contamination rates 

In the Gaia processing pipeline a source will only be charac- 
terised as an ultra-cool dwarf if the following conditions are met: 

1 . Its estimated distance is smaller than the maximum distance 
at which an ultra-cool star can be detected. This is the dis- 
tance at which the brightest ultra-cool star would have a Gaia 
limiting magnitude of G=20. We estimated the maximum 
brightness of an ultra-cool star from the BT-Settl models, 
and it corresponds to the hottest model with the lowest sur- 
face gravity. 

2. The source is fainter (in the G band) than the biightests 
model placed at the same distance as estimated by Gaia. 

3. The BP-RP colour index is higher (redder) than the mini- 
mum (bluest) colour index found in the model libraries. 

4. The celestial coordinates and proper motions are not consis- 
tent with solar system Keplerian motions. 
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Table 1: Predicted number of counts assuming solar metallicy and auxiliary relations. 

Coun tSchKOi CoMilfiBT-Seltl T'cff (Gbp-Grp) 



spectral type Mj 



M3-4 V 


9.07 


10.24 


1.6 X 10' 


5.7 X 10'' 


3269 


2.86 


M4-5 V 


10.14 


10.86 


5.4 X lO*" 


2.4 X 10* 


3065 


3.33 


M5-6 V 


11.91 


11.39 


3.4 X 10' 


7.0 X 10' 


2889 


3.74 


M6-7 V 


12.90 


11.92 


4.3 X 10" 


1.7 X 10' 


2731 


4.10 


M7-8 V 


14.16 


12.47 


3.0 X 10^ 


3.0 X 10" 


2583 


4.44 


M8-9 V 


14.68 


13.05 


1.0 X 10^ 


9.4 X 10^ 


2441 


4.77 


M9-L0V 


15.23 


13.66 


3.0 X 10^ 


2.5 X 10^ 


2303 


5.03 


LO-1 V 


15.27 


14.30 


2.3 X 10^ 


8.4 X 10^ 


2166 


5.07 


L1-2V 


15.69 


14.95 


1.7 X 10^ 


4.9 X 10^ 


2032 


5.10 


L2-3V 


16.19 


15.44 


1.0 X 10^ 


2.9 X 10^ 


1902 


5.62 


L3-4V 


16.67 


15.87 


5.4 X 10' 


1.3 X 10^ 


1777 


5.35 


L4-5V 


17.14 


16.38 


3.2 X 10' 


7.2 X 10' 


1660 


5.15 


L5-6V 


17.61 


16.97 


1.8 X 10' 


4.3 X 10' 


1554 


8.02 


L6-7V 


18.07 


17.57 


8.7 


1.8 X 10' 


1460 


9.05 


L7-8V 


18.52 


18.09 


3.7 


6.8 


1381 


9.85 


L8-9V 


18.90 


18.49 


1.5 


2.7 


1318 


10.38 


L9-T0 V 


18.95 


18.78 


9.4 X 10"' 


1.1 


1270 


10.81 


TO-1 V 


19.01 


18.97 


6.1 X 10"' 


6.1 X 10"' 


1237 


11.17 


T1-2V 


19.04 


19.10 


6.4 X 10"' 


5.5 X 10"' 


1214 


11.44 


T2-3V 


19.04 


19.19 


1.1 


8.0 X 10"' 


1196 


11.64 


T3-4V 


19.07 


19.29 


2.0 


1.3 


1176 


11.84 


T4-5V 


19.22 


19.44 


2.9 


1.9 


1144 


12.15 


T5-6V 


19.61 


19.68 


2.6 


2.1 


1086 


12.61 


16-1 W 


20.35 


20.05 


1.2 


1.8 


986 


12.99 


T7-8V 


21.60 


20.62 


4.0 X 10"' 


1.4 


824 


12.46 
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Fig. 4: Predicted number of counts per spectral type bin (a) and per apparent G magnitude (b), in logarithmic units. The black line 
in the left panel coiTesponds to the derivation based on the relation between spectral type and I-band absolute magnitude included 
in Caballero et al. (2008), while the blue continuous line corresponds to the relation derived from the BT-Settl model family and 
the SLC calibration. The two horizontal (dashed) lines indicate the levels of predicted counts equal to one and ten. The right-hand 
side plot has been obtained assuming the relation between effective temperature and I-band absolute magnitude derived from the 
BT-Settl models and the SLC calibration. 



Sources that fulfil these criteria (within some margin that de- 
pends on the measurement uncertainties) are subsequently pro- 
cessed to estimate their effective temperature and surface gravity 
as described in section 3. Objects not detected in the BP band 
(but fulfilling all other criteria) will nevertheless be selected be- 
cause the non-detection is itself an indication of a red spectrum. 
This could imply a potential contamination of the UCD sample 
by faint, nearby objects close to the Gaia limiting magnitude and 
colour indices that are positive (thus excluding white dwarfs), 
but bluer than the bluest UCD. The software that implements the 



analysis of ultra-cool dwarfs distinguishes between the selection 
module that decides whether to process a given Gaia source, and 
the processing module that estimates the effective temperature 
and gravity of the selected sources. 

To define the boundaries of brightness, colour index, and dis- 
tance that define the region of ultra-cool sources, we analysed the 
GOG simulations of the AMES-Cond, AMES-Dusty and BT- 
Settl model libraries (although the latter is prefered, and the for- 
mer have not been implemented in the selection software). Using 
GOG, we found that the most distant ultra-cool object in the li- 
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Fig. 3: Maximum distances at which an ultra-cool dwarf can be 
detected by Gaia at the limiting magnitud G=20 as a function 
of its absolute magnitude in the I band. These have been derived 
from BT-Settl models (filled circles) and the continuous lines 
represent the interpolation used in deriving the expected counts 
per spectral type bin in Table 1 . The black continuous line cor- 
responds to log(g)=5.0 and the blue line to log(g)=3.5. The top 
axis shows the effective temperature measured in Kelvin for a 
log(g)=5.0 object with the absolute I magnitudes shown in the x 
axis, according to the BT-Settl models. The Tf-g-Mi mapping is 
only bi-valued below 600 K. 



brary (a source with Teff - 2500 K and log(^) = 3.5 at G=20) cor- 
responds to a distance of 373 pc. Thus, the first criterion limits 
the processing of sources by the UCD module to sources within 
this radius (this is necessary optimise the processing time per 
source). The second criterion examines the brightness of each 
detected source and compares it with the aforementioned bright- 
est model placed at the same distance. We expect a strong con- 
tamination from main-sequence stars with temperatures above 
the 2500 K limit because low-gravity sources are brighter than 
main-sequence ones for temperatures above 1600 K. In fact, we 
derived from the GOG simulation of the model libraries that a 
dwarf star of 4180 K has the same G magnitude as the lowest 
gravity, hottest ultra-cool giant considered in this work (Teff - 
2500 K, \og{g) = 3.5). Thus, main-sequence objects up to this 
temperature will be selected by the selection module according 
to this criterium despite their higher temperatures (although the 
colour index criterium will reject a fraction of these hotter stars; 
see below). We could have used the thresholds corresponding 
to the main-sequence objects, but this would have resulted in 
the potential loss of low-gravity sources in the selection process. 
We opted for an inclusive set of criteria in spite of the expected 
high contamination rates by dwarfs hotter than 2500 K. To es- 
timate the contamination rate (caused by the criteria but also 
to the measurement errors), we conducted a numerical experi- 
ment in which we populated a one-kiloparsec cube with sources 
uniformly distributed in space (this is much larger than actually 
needed in view of the estimated maximum distance to an UCD). 
The total number of sources in each temperature bin was gener- 
ated using a probability density function derived by interpolat- 
ing and normalising the volume densitites tabulated in Caballero 
et al. (2008). Since these correspond to main-sequence sources. 



we neglected the contribution from ultra-cool giants in the fol- 
lowing estimates. 

For each star, and given the effective temperature and the 
distance to the centre of the cube, we computed the apparent 
G magnitude and the (V-I) colour index by interpolating in the 
BT-Settl library restricted to dwarfs with values of log(§) x 5.5. 
We subsequently used the current estimates of the uncertainties 
in the Gaia magnitudes and parallaxes (see Jordi et al. (2010) 
and de Bruijne (2009) respectively) to generate mock measure- 
ments of the distances and G magnitudes by sampling from nor- 
mal distributions with the prescribed uncertainties. We did not 
use a full covariance matrix since none was available at the time 
of writing. We show in Fig. 5 the error models that we used in 
the simulations. The uncertainty estimate in the G magnitude is 
based on Eq. 6 in Jordi et al. (2010) and on the mission param- 
eters available at the time of writing. We included a calibration 
error cTcai =30 mmag for a single transit. This value is only an 
educated guess since the final cTcai can only be estimated during 
the operational phase. The uncertainty in the parallax is based 
on Eq. 1 in de Bruijne (2009) (see de Bruijne (2012) for a more 
recent review without analytic expressions for the astrometric 
uncertainty). 




o 

« 
o 



CJI 

o 



Fig. 5: Current estimates of the end-of-mission uncertainties in 
the measurements of the G apparent magnitude (black line) and 
the parallax (dashed and dotted lines) as a function of the G ap- 
parent magnitude. The dashed line corresponds to (V-I) - 4 and 
the dotted line to (V-I) = 7.5, a plausible range for the V-I colour 
index according to the model libraries. 



The resulting sample of stars in the kiloparsec cube was 
then examined to determine the properties of the subsample 
that fulfilled the criteria enumerated above (except for the non- 
Keplerian motion). We found that none of the ultra-cool dwarfs 
thus generated is missed in the selection process due to errors in 
the measurements. This was expected given i) the inclusiveness 
of the criteria, ii) the relative proximity of these faint sources to 
the Sun, and iii) the fact that we only generated main-sequence 
objects in the simulation (and not low-gravity young sources that 
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would be closer to the selection thresholds; this question will be 
re-examined during the software validation phase with real Gaia 
data and may result in more inclusive thresholds if we identify 
examples of this kind missed by the selection module). The un- 
certainties in the G magnitude and colour index typical of UCDs 
according to the simulations are negligible in this context. 

The (Gbp-Grp) colour index thresholds in the selection crite- 
ria are defined by the bluest model in the GOG simulations of the 
model libraries. For the BT-Settl library of models and metallici- 
ties [M/H] between -2 and 0.5, we find that the (Gbp-Grp) colour 
index is in the 4.1-14.5 range. The distribution of effective tem- 
peratures in the resulting sample of non-UCD stars that fulfil the 
criteria is shown in Fig. 6. The shape is determined by the com- 
bination of the various selection criteria and the uncertainties in 
the measurements. 



IIIININININ^ 

2500 2600 2700 2800 2900 



Fig. 6: Histogram of the temperatures of non-UCD stars that pass 
the selection criteria of the Gaia UCD module. 



Fig. 6 impUes a contamination rate of abound 85%. This is 
explained by the exponential increase in volume densities for 
the mid- to late-M spectral types. To take this high contamina- 
tion rate into account, we generated regression models in Sect. 3 
that are capable of predicting efi'ective temperatures up to 4000 
K, and not only in the UCD regime below 2500 K. Therefore, 
even if the selection module selects these contaminants, the pro- 
cessing module will allow for the obtention of purer samples of 
UCDs by filtering out objects with estimated temperatures above 
a given threshold. 

In computing the contamination rate we did not consider 
extra-galactic contaminants because their parallax measure- 
ments would be incompatible with the selection criteria defined 
above. 

3. Methodology 

As mentioned in Sect. 1, Gaia is expected to detect and process a 
number of sources close to one billion. This order of magnitude 
necessitates automating the detection and characterization pro- 
cesses. The main objective of the eighth coordination unit of the 
Gaia DPAC is to implement and evaluate automated procedmes 
to derive the astrophysical parameters of the sources detected 
by the Gaia instruments. This necessarily involves techniques 
developed in the fields of statistical learning and data mining. 
In the framework of statistical learning methods, this is accom- 
pUshed with regression models that are constructed from sets of 



examples (the so-called training set) that link the independent 
variables (the Gaia observations in our case) with the dependent 
variables that we wish to infer. In this section we describe sev- 
eral techniques taken from the field of statistical learning, aimed 
at providing a reproducible, systematic characterization of the 
Gaia UCD candidate sources. 

In collecting the aforementioned sets of examples, libraries 
of stellar models and associated synthetic spectra provide a ho- 
mogeneous and consistent set of examples that uniformly cover 
the parameter spaces under consideration. In contrast, catalogues 
of observed spectra with associated astrophysical parameters de- 
rived from them tend to be fragmentary in nature: each collection 
covers only a reduced range of parameters, and it is necessary 
to combine several catalogues to obtain a coverage that, even 
then, may contain gaps and/or insufliciently sampled regions 
(this is especially true for the parameter log(^), which is miss- 
ing from most compilations). Each catalogue has its own ob- 
servational setup, systematic errors, and selection effects. Also, 
observational biases may favour the abundance of examples in 
particular regions of the parameter space, which then translate 
into systematic biases in the predictions. Furthermore, parame- 
ter estimation across catalogues can be inconsistent, and we in- 
deed find slightly difi'erent spectral types assigned to the same 
source in different catalogues, which reveals a certain degree of 
subjectivity in the assignements. With all these considerations 
in mind (reproducibility, consistency, homogeneity, and unifor- 
mity of the training set), we prefered to construct our regres- 
sion models from the aforementioned synthetic libraries. These 
libraries parameterise the models with physical magnitudes (ef- 
fective temperatures, gravities, and metalUcities) and not with 
spectral types. 

Model libraries are nevertheless imperfect in the sense that 
they do not reproduce each and every spectral feature in the real 
spectra of UCDs or its exact dependence on the physical pa- 
rameters. The problem of the potential mismatch between stellar 
model libraries and observations appears whenever physical pa- 
rameters are to be produced, beyond phenomenological descrip- 
tions such as spectral types. Spectral types can be inferred with- 
out the need for synthetic models, while going from spectra (or 
spectral types) to physical parameters requires model libraries to 
interpret them. In this work we do not attempt to build regression 
models to infer spectral types. Given the low spectral resolution 
of the Gaia red spectrophotometer (RP), and its spectral cov- 
erage, most of the spectral features used to decide the spectral 
type remain unresolved or unobserved, and thus, spectral types 
derived from them would be of limited use. 

The systematic efl'ects of the selection of model libraries can 
be characterised to some extent by comparing the predictions 
from different model families like those presented in the previ- 
ous section (AMES-Cond, AMES-Dusty, and BT-settl). In Sect. 
4 we also attempt to validate our models with an external set of 
effective temperatures derived from ground-based observations 
and spectral types via a calibration that, inevitably, encompasses 
another synthetic stellar library. The two families of models in- 
troduced in the previous section (AMES-Cond and -Dusty, and 
the BT-Settl models) are used to define the relationship between 
Gaia observables (the RP spectrum) and the parameters that we 
intend to estimate, namely T^ff and log(g). This relationship is 
captured in a regression model (not to be mistaken for the physi- 
cal models of the stellar atmospheres and the resulting synthetic 
spectra) constructed from the set of examples defined by the two 
synthetic spectral libraries. Each spectrum in the libraries (e.g. 
those represented in the lower row of Fig. 1) is characterised 
by the set of physical parameters that identifies the stellar at- 
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mosphere used to synthesize it. This set of examples composed 
of the spectrum plus the corresponding physical parameters is 
referred to in the following as the training set. The set of inde- 
pendent examples used to assess the accuracy of the models will 
be referred to as the test set. 

As mentioned before, the training set is constructed using the 
AMES-Cond, -Dusty, and BT-Settl model libraries (restricted to 
effective temperatures below 4000 K)) and simulating Gaia ob- 
servations of these synthetic models using GOG (version 8.0). 
Gaia has a nonninal duration of five years. During this period, 
a source is observed on average approximately 70 times. Each 
time a source is observed (i.e., each transit), an epoch RP spec- 
trum is obtained. The characteristics of this spectrum depend on 
the instrument design including prism and CCDs, and the transit 
geometric properties. In general, we can assume that the spec- 
trum will be spread along 60 spectral bands with a non-uniform 
dispersion. In each transit, the position of the source (continu- 
ous) spectrum may vary with respect to the discrete CCD pixel 
array depending on the transit geometric details. This transit de- 
pendence results in diff'erent spectra obtained for each epoch, 
because the wavelength coverage of each pixel is different. This 
subresolution information can then be used to produce oversam- 
pled combined spectra. In Fig. 1 and in the experiments carried 
out in this work that are described in the following sections we 
assumed an end-of-mission oversampling factor of three, result- 
ing in a total number of flux bins of 180. 

3.1. Regression models 

In this section we describe three different types of statisti- 
cal regression models: A:-nearest neighbours (Cover & Hart 
1967), Gaussian processes (see e.g. Bishop 2006; Rasmussen 
& Williams 2006), and Bayesian inference (see e.g. Sivia & 
Skilhng 2006). A previous analysis also included support vec- 
tor machines (Vapnik 1995; Cortes & Vapnik 1995) and kernel 
partial least squares (KPLS, Rosipal et al. 2001). 

A:-nearest neighbours estimation (kNN) is by far the sim- 
plest model and derives the parameter values as the weighted 
average of the elements in the training set that are closest to the 
input spectrum in a given metric. In our case, the euclidean dis- 
tance is used to define proximity and the weights are defined as 
the inverse of this distance. A maximum-hkehhood estimate of 
the uncertainties in the estimated parameters can be computed 
by modelling the distribution of inverse distances to the nearest 
neighbours (under the assumptions that the model grid is suf- 
ficiently dense and the neighbourhood has the appropriate size 
to sample a unimodal likelihood). It has the disadvantage that 
the full training set has to be stored and accessed each time that 
the regression model is used to predict the physical parameters 
iTgg and log(g)) of a source, and furthermore, it is severely af- 
fected by the so-called curse of dimensionaUty. This is reflected 
in the need for exponentially growing training set sizes as the 
problem input dimensionality increases. The exponential growth 
is required to ensure that the nearest neighbours are sufficiently 
close to provide an accurate estimate of the parameters. 

Support vector machines (SVMs) and Gaussian processess 
(GPs) are two examples of kernel methods. These methods trans- 
form the regression problem into a dual representation where the 
constituents of the model are no longer the input features (or 
nonlinear mappings thereof) but their scalar products expressed 
as kernel functions. More details of this dual representation can 
be found in the textbooks by Bishop (2006) and Hastie et al. 
(2001). 



A GP is defined as a probability distribution over functions 
such that the joint probabihty of the random variables defined 
by their evaluations at a certain set of input vectors (the train- 
ing set in our case) is Gaussian. If we assume that the nature of 
our problem is such that the probabilistic distribution for both 
the targets of our training set and the Gaia observations to be 
characterised is well captured by a multivariate Gaussian distri- 
bution, we can construct the model by computing its covariance 
matrix. It turns out that we can calculate it by evaluating the ker- 
nel functions at the input vectors of the training set. These kernel 
functions encode both the assumed error model for our determi- 
nation of the target values in the training set and the length scale 
for the correlations between the examples in the training set. In 
the full Bayesian treatment these two parameters can be consid- 
ered hyperparameters and are marginalised out. Here we have 
determined optimal values for them through exhaustive cross- 
validation experiments (see below). For the SVM, the model 
representation is sparse in the sense that it only depends on ker- 
nel evaluations on a small subset of the training set (the support 
vectors). In regression problems, the support vectors are those 
that he within the boundaries of, or outside a so-called insensi- 
tive tube. Support vector regression involves the search for best 
(minimum error) solution in a space of hyperparameters, similar 
in nature to those discussed for GPs. 

In this section we describe a paraUel aproach to parameter 
estimation based on Bayes' theorem. The main advantage of 
Bayesian parameter estimation stems from the fact that it pro- 
vides not only an estimate of the parameters but also a full mul- 
tivariate probability density distribution for the set of parameters 
given the observations. 

We assume that we observe the spectrum s of an ultra-cool 
dwarf. In the Bayesian framework, we seek the probability den- 
sity function (PDF) of the physical parameters of the UCD, given 
this spectrum s. Bayes' theorem provides this PDF as 



nie\s) 



mis) 



(1) 



where 9 - (Tg^, log(^)) is the vector of parameters that we 
intend to derive, niO) is the prior probability distribution of these 
parameters, and feis) is the UkeUhood of the spectrum given the 
parameters. 

The denominator in the right-hand side of Eq. 1 m(s) is the 
prior predictive distribution or evidence defined as 



«(s)= fe(s)n(e)de. 



(2) 



and can be viewed as a normalisation constant. 

The likelihood term encompasses a predictive model for the 
spectra given the parameters, together with a probabilistic er- 
ror model. In our case, we used the two stellar hbraries men- 
tioned above (COND-i-DUSTY and BT-Settl) to build the pre- 
dictive model. This is captured in two three-layer perceptrons 
(a kind of neural network) each trained with the stellar spec- 
tra of the corresponding Ubrary. The neural network captures a 
multivariate regression model that can be used to perform an in- 
terpolation. Extensive experiments to derive the optimal network 
architecture result in hidden layers of 20 hidden units for the two 
model libraries. With it, we can generate output synthetic spectra 
for any input value of 6 within the grid boundaries. The neural 
network exactly predicts the spectra in the COND-i-DUSTY and 
BT-Settl grids, and smoothly interpolates for intermediate val- 
ues. The unique flux predicted by the neural network for each 
wavelength in the spectrum can be turned into a probabilistic 
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statement by convolving it with the predicted measurement er- 
rors in the current Gaia model. At present, the error model for 
the RP spectra consists of a Gaussian distribution with a covari- 
ance 2 that depends on a series of instrument parameters, on the 
number of transits, and on the flux itself (see Jordi et al. (2010), 
discussed above). Thus, 

/e(s) = /(x„E)(s), (3) 

with Xg being the prediction from the neural network for the pa- 
rameter set 6, and /(xj,,i;)(s), the probability density function of a 
normal distribution Nixg, S) evaluated at the observed spectriun 
s. 

Nested sampling (Skilling 2006) is a Monte Carlo procedure 
used to calculate m(s) that can extract samples from n{6\s) as a 
by-product of these calculations. The algorithm exploits the rela- 
tion between the UkeUhood /e(s) and the prior volume X defined 
by 

XiA)= f n{d)de, (4) 

J/e(s)>^ 

i.e., the volume of n{d) over the region of parameter space con- 
tained within the A likelihood iso-contour /©(s) = A. Our prior 
density n{-) can be set in two ways: i) according to a bidimen- 
sional uniform distribution in (400, 3500) X (3.5, 5.5), or ii) using 
a physical prior based on the temperature histograms discussed 
in the previous section, and a log(^) prior that favours values typ- 
ical for the field main-sequence dwarfs. We here studied the pre- 
dictions of a Bayesian module for regression based on a flat prior 
in both parameters Teff and \og(g). This represents the frame- 
work for the analysis of the influence of physical priors (based 
for example on estimated detection rates such as those presented 
in section 2 or on the astrometric and photometric data provided 
by Gaia) to be included in a subsequent paper in this series. 

We used ellipsoidal sampling (Shaw et al. 2007) to estimate 
the parameters Teff and log(^). This method is a variant of nested 
sampling that approximates the iso-hkehhood of the point to be 
replaced in the nested sampling step by an n-dimensional el- 
lipsoid (e/v) derived from the covariance matrix of the current 
active points (see Skilling (2006) for more details). Since our 
stellar models are restricted in parameter space to the region 
500 < Teff < 3500, 3.5 < log^ < 5.5, points drawn from 
but outside of the parameter space will be discarded. 

Even though ellipsoidal sampling is aimed at simplifying the 
computation of the evidence, it is also possible to derive poste- 
rior probabilities from it as a by-product. Once the algorithm has 
converged, the resulting sample can be interpreted as a sample 
from the posterior probability if we weight the importance of 
each point by a factor defined as 



with w,- = l/2(X',_i - X,+i). The justification for this weight- 
ing scheme is described in Skilling (2006). Thereafter, we can 
obtain summary statistics such as the mean of the posterior prob- 
abihty density of eff'ective temperatures using the classical first- 
order moment of the distribution: 

n 

TeB = ^ Teffj ■ Pi, 

!=1 

with Teif,! the effective temperature value in the ?-th sample, and 
Pi, the weights calculated according to Eq. 5. We use this sum- 
mary statistic (the first-order moment) to estimate T^g, and will 



discuss the effect of using other alternatives such as the posterior 
mode in a subsequent paper. It is planned that the Gaia catalogue 
will contain the samples obtained for each UCD candidate. 

Ellipsoidal sampling was preferred over simpler algorithms 
such as Metropolis-Hastings or nested samphng because tests 
carried out with these algorithms and several proposal densi- 
ties produced too low acceptance rates because of the particular 
shape of likelihood landscape. This results in posterior samples 
that are not independent. These problems are more severe in on- 
going appUcations of the software to highly multi-dimensional 
related problems such as the estimation of star-i-disk parameters 
in pre-main-sequence systems. 

3.2. Preprocessing 

Before constructing the models, the spectra were normalised to 
yield an area equal to 1 to have an appropriate scale of values 
that is robust to noise and isolated outUers. This normalisation 
removes the information relative to the integrated energy flux 
that, combined with the distance measured by Gaia, can pro- 
vide indications of the Jeff and log(g) values of the source. This 
information was incorporated at a later stage, and also in the 
consistency checks. 

Other preprocessing steps were explored to determine their 
impact on the performance of the algorithms. In particular, we 
tested denoisification strategies based on wavelet decomposition 
and moving averages. The GOG allows generating noisy simu- 
lated spectra for any number of transits and apparent magnitude. 
In the experimental setup we included GOG simulations for a set 
of four apparent G magnitudes (G=15, 18, 19 and 20) and two 
number of transits (28 and 70). Seventy transits is an estimate 
of the average number of transits after five years of observations 
(the nominal Gaia lifetime), whilst 28 transits corresponds to the 
average after two years of observations, when a reliable evalua- 
tion of the algorithms can be attempted. The number of transits 
strongly depends on the position of the source on the sky (see 
e.g. Lindegren et al. 2012). 

Several wavelet bases were tried in the experiments includ- 
ing several orders of the Daubechies, Coiflet, best-located and 
least asymmetric wavelets. In addition to denoisification strate- 
gies, we explored two data compression approaches: the well- 
known principal component analysis (Pearson 1901) and the lo- 
cal approach based on difliision maps (Coifman & Lafon 2006). 

3.3. Internal validation of the regression models 

The experiments described in this section were carried out with 
regression models trained with the so-called nominal dataset. 
This comprised GOG simulations of the original spectra in the 
synthetic spectral Ubraries (restricted to efi'ective temperatures 
below 4000 K), and it corresponds to the nodes of a grid in 
the space of Tg^ and \og(g), with grid spacings of 100 K in Tetf 
and 0.5 dex in log(g). The performance of the regression mod- 
els was measured by analysing the distribution of residual errors 
in the so-caUed random dataset. This comprises spectra linearly 
interpolated from the nominal grid at values of the physical pa- 
rameters (Tgff and log(g)) randomly spread in the evolutionary 
tracks provided by the authors of the libraries of stellar models. 
The regression models can be understood as continuous nonlin- 
ear mappings between the 180-dimensional space of spectra and 
the one-dimensional space of the parameter under consideration. 
The random dataset is not an independent test set because it is 
derived from the nominal dataset using multi-linear interpola- 
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tion. In this sense, the peiformance evaluation measures not only 
the abiUty of the regression models to reproduce the training set, 
but unfortunately also the fidelity to the multi-linear interpola- 
tion between the grid nodes. We refer to the performance mea- 
sures described in this section as intemal errors, in the sense that 
they measure the ability of the regression models to reproduce 
the mapping between spectra and astrophysical parameters that 
is inherent to the training set in the UCD domain (even if the per- 
formance also measures the fidelity to the multi-linear interpo- 
lation). Internal errors, thus, evaluate only the regression model 
and its robustness against noise, but not the validity of the train- 
ing sets or their ability to reproduce observed spectral features. 
In Sect. 4 we analyse the performance of the regression mod- 
els by applying them to observed spectra (downgraded to the 
spectral resolution of the Gaia instruments). These performance 
measures, thus, comprise the intemal errors and the ability of the 
triiining sets to reproduce the spectra of real UCDs. 

The complete description and analysis of the preliminary set 
of experiments described in the previous paragraphs is beyond 
the scope of this article, but the conclusions derived from it can 
be summarised as follows: 

- Residual errors increase near the boundaries of the training 
sets. The region below 500 K was particularly problematic 
because prediction errors were unacceptably large and in- 
cluding it in the training sets also degraded the performance 
for higher temperatures. It was therefore removed from the 
training sets. 

- Denoisification is best achieved by using a moving average 
filter as measured by the root mean squared reconstruction 
error (RMSE), and the mean and median reconstruction er- 
ror. 

- The lowest RMSE obtained with wavelet denoising is 
achieved with the Coiflet mother wavelet of the order of 18; 
similar RMSE are obtained with the best located wavelet of 
the order of 14. 

- Training the algorithms with a noisy training set adjusted to 
the signal-to-noise ratios of the test spectrum yields better 
results than training with noiseless spectra and denoisifying 
the test spectrum, regardless of the denoisification strategy 
(moving averages or wavelets). 

- Prediction errors obtained with KPLS models are signifi- 
cantly larger (RMSE systematically above 200 K for all G 
magnitudes and numbers of transits) than those obtained 
with GP, SVMs, or kNN. models; 

- Prediction errors from cross- vaUdation experiments with the 
various models (except KPLS) and noiseless spectra are all 
the same within the experimental uncertainty as measured by 
the standard deviation of the RMSE sample obtained from 
ten experiments of ten-fold cross validation. The RMSE er- 
ror for the noiseless training set is 7 K. 

- kNN models systematically outperform both SVM and GP 
models when applied to noisy spectra in the estimation of 
Teff values. Figure 7 shows a typical case where both the 
training (nominal) set and the test (random) set were sim- 
ulated with noise properties for G=20 and 70 transits. The 
left column shows the distribution of residuals for the GP 
model (reff and log(^) in the upper and lower panel, respec- 
tively) and the right column shows the same scatter plots for 
the kNN model. The RMSE in the T^fs predictions is 10 K 
for the kNN model and 62.3 for the GP model. The RMSE 
of the kNN model is unrealistically low as shown in the next 
section, and due in part to the high density of examples in the 
training set, in this context of cross-validation experiments to 



derive internal errors. In this setup, where we assumed that 
the training set of synthetic spectra reproduces the expected 
observed spectra well , the nearest neighbours are as close as 
allowed by the density of the grid of training examples (given 
the relatively high signal-to-noise ratio of the end-of-mission 
spectra). The performance of the SVM model is overall very 
similar to that of the GP model. 
- Reducing the input space dimensionality with PCA (pre- 
serving 95% of the variance) deteriorates the model perfor- 
mances by up to 100% with respect to the complete input 
spectra. This is so except for the lower signal-to-noise ra- 
tios (G=20) and in the combination of PCA and kNN, where 
the RMSE improves by a 30%. The improvement decreases 
as the noise diminishes and for G=15, the two models (with 
and without PCA compression) show the same RMSE. Data 
compression with the nonlinear technique known as diffu- 
sion maps (Coifman & Lafon 2006) results in similar perfor- 
mances as the PCA approach at a much higher computational 
cost. 

In Sect. 4 we concentrate on the external validation of the 
kNN and GP models only, because the SVM model perfor- 
mances are remarkably similar to the GP models, but do not 
provide estimates of the prediction uncertainty in a straightfor- 
ward way. As a result of the previous considerations, we have 
several regression models based on the A:-nearest neighbours al- 
gorithm and GP. There is a version of each model for the 28 
and 70 transits cases, and for the various signal-to-noise ratios 
that correspond to values of the G magnitude G=15, 18, 19, and 
20. For each of these cases we have a model trained in the full 
input space, or on a PCA compressed version of it. All param- 
eters of each model, such as the number of nearest neighbours 
or the kernel and noise parameters for the GP models, are deter- 
mined using ten-fold cross-validation of the models obtained in 
an exhaustive exploration of the parameter space. The resulting 
models are used in section 4 to analyse the expected uncertain- 
ties when applied to spectra of ultra-cool dwarfs observed from 
ground and simulated with GOG. 

Unfortunately, the intemal validation described above for the 
kNN and GP models overestimates the real accuracy because it 
does not include several systematic eff"ects. The physical models 
used in constructing of the training set (AMES-Cond, AMES- 
Dusty, and BT-Settl) do not exactly reproduce all spectral fea- 
tures and their correlation with the physical parameters encoun- 
tered in reality. This may be due to several reasons, such as in- 
complete line lists, inaccurate Une/band opacities, or mathemat- 
ical simplifications introduced for the sake of tractability (e.g. 
in the dust cloud formation and convection, in the diffusion of 
chemical species, in the departures from equilibrium between 
the dust and gas phases, or in the neglect of rotation). To ac- 
count for these systematic errors, a second battery of experi- 
ments is discussed in section 4 where ground-based spectra of 
well-known ultra-cool stars are degraded to the Gaia resolution 
and are convolved with its spectral response and error model us- 
ing GOG. We aim at estimating the error that affects each of the 
regression models in the various T^s and log(^) ranges covered 
by the empirical spectral libraries. The spectra used for the ex- 
ternal validation of the regression models are described in detail 
in the next section. 

4. Validation with real spectra 

In this section, we use ground-based spectra to estimate the to- 
tal errors of the regression models. These total errors include 
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Fig. 7: Cross-validation errors for end-of-mission spectra of UCD stars at G-20. The upper panels show the eiTors in the T^ff 
estimates for the GP model (left) and the ^-nearest neighbours model (right). The lower panels show the corresponding errors in 
the \og(g) estimates. The colour code for the error scale is shown in the right-hand side of each row. Effective temperatures are 
measured in Kelvin and gravities in cm/s^. 



the internal validation errors described in the previous section, 
and errors due to the imperfect representation of real spectra 
by synthetic models. The ground-based spectra used for vali- 
dating the regression models are the compilations by A. Reid^ 
(Burgasser et al. 2000; Delfosse et al. 1997; Gizis et al. 2000b,a; 
Kirkpatrick et al. 1999, 2000; Reid et al. 1999, 2000; Strauss 
et al. 1999; Tsvetanov et al. 2000), S. Leggett^ (Chiu et al. 
2006; Golimowski et al. 2004; Knapp et al. 2004), the NIRSPEC 
Brown Dwarf Spectroscopic Survey"* (McLean et al. 2003), and 
the IRTF spectral library^ (Gushing et al. 2005; Rayner et al. 
2009). None of the observed spectra covers the full Gaia wave- 
length range, which is especially due to the lack of optical ob- 
sei-vations at the coolest end. Therefore, it was necessary to 
complete them with models before simulating the Gaia obser- 
vations of these stars for ad hoc G magnitudes (and therefore 
noise levels) and total number of transits. This was accomplished 
by degrading the resolution of the model spectra to that of the 
observed spectrum and finding the model that yields the mini- 
mum fit in the wavelength range of overlap restricted to the 
RP passband (after removing artefacts and wavelength ranges 
with poor signal-to-noise ratios). This model was then used to 
complete the observed spectrum. We used all available models 
from the COND, DUSTY, and BT-Settl libraries to find the best- 
fitting model to each of the observations, and in all cases the BT- 
Settl model produced ax^ statistic superior or comparable to the 
COND/DUSTY models. Spectra with gaps, small overlap with 
the RP band, or very low signal-to-noise ratios were removed 
from this empirical validation set to avoid biases. 

^ http : / /www . stsci . edu/ ~ inr /ul tracool . html 
^ http : //staf f . gemini . edu/- sleggett/LTdata .html 
* http : / /www . astro . ucla . edu/~mclean/BDSSarchive/ 
#lowres 

' http://irt£web . ifa.hawaii . edu/-spex/IRTF_Spectral_ 
Library/ index . html 



The best x^ fits to the full resolution spectra give us a lower 
limit to the eiTors attainable in the Gaia Teff estimates. Figure 8 
shows in the y-axis the effective temperatures derived from the 
best x^ fits described in the previous paragraph as a function 
of the spectral types assigned in the libraries. We introduced a 
small jitter (characterised by a standard deviation equal to 30K) 
in the values of the fitted temperature to enhance the visibility of 
the stars with the same spectral type. In the same plot we have 
included the calibration of effective temperatures with spectral 
types by Stephens et al. (2009) as a continuous line. It corre- 
sponds to the calibration derived from optical spectral types of 
L dwarfs and infrared spectral types of the T dwarfs, valid in the 
M6 to T8 range (the SLC calibration). The SLC calibration relies 
on model atmospheres described in Saumon & Marley (2008), 
and references therein. The dashed lines correspond to the same 
calibration plus/minus 250 K. The right-hand panel shows the 
difference between the predicted temperature and the tempera- 
ture obtained applying the SLC calibration to the spectral types 
assigned in the compilations, as a function of the latter. The 
scatter of T^ff values in Fig. 8 around the SLC calibration re- 
flects errors of the regression process, but also the uncertainty 
in the calibration between spectral types and effective tempera- 
tures. In this work, we are taking the effective temperature de- 
rived through the SLC relation as the ground truth. However, 
the SLC relation, derived from high-resolution infrared spectra, 
which are much more appropriate for determining an effective 
temperature than the Gaia spectrophotometry, has an intrinsic 
scatter estimated by the authors as IQQK (although the article 
lacks details about how the scatter was estimated, and visual in- 
spection of the plots suggests a much higher scatter at least in the 
L sequence). Furthermore, Reyle et al. (2011) compared a cali- 
bration derived from the BT-Settl model library (the one used in 
this work) with several other calibrations for the M spectral sub- 
sequence in the literature, and found variations of the order of 
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200-300 K for each spectral subtype. Bayo et al. (201 1) showed 
even higher values of the scatter for the same spectral type range 
(they derived the T^g values from global fits to the SED). This 
shows that the calibration of the relationship between effective 
temperatures and spectral types is a problem that is not yet fully 
solved, and that the exact Teg values derived for an UCD de- 
pend on the calibration used, with typical uncertainties of a few 
hundred Kelvin. The predictions of the regression modules pre- 
sented in this section have to be judged in this context. In this 
respect, our regression models for the Gaia data have the ad- 
vantage that they will provide UCD temperatures derived con- 
sistently from the BT-Settl model family. 

Three stars lacked precise spectral subtypes in the 
spectral compilations, namely 2MASS 1237392+652615, 
SDSS 1346464-003 150, and SDSS1624144-H002916. We as- 
sumed them to be T6.5 (Kirkpatrick et al. 2011), T6.5 (Reyle 
et al. 2010), and T6 (Kirkpatrick et al. 2011). Figure 8 shows 
a tendency to assign effective temperatures around 1700 K for 
sources with spectral type L. 

Table 2 lists the mean difference // between the T^g esti- 
mates and the effective temperatures derived from the SLC cal- 
ibration and the spectral types cited in the spectral compilations 
(hereafter bias); the standard deviation with respect to the bias- 
corrected mean (cr) displayed by the four compilations; and the 
RMSE without correcting for the mean bias. We list these values 
for the two model famiUes (AMES Cond + Dusty and BT-Settl) 
used in this work. All RMSE values are given in Kelvin. 

The compilations by Reid and Leggett have a good overall 
wavelength coverage in the RP range, and spectra with poor cov- 
erage were removed from the validation set. Since the expected 
emission in the wavelength regions of the RP range not covered 
by the observed spectra (bluewards of 750 nm) is negligible in 
the Teff and log(g) parameter ranges under consideration, we ex- 
pect the completions to have Uttle or no effect on the subsequent 
parameter estimation with the models described in the previous 
section. The NIRSPEC library is in the opposite case, with a 
majority of spectra starting around 1.1 /dm. At this wavelength, 
the RP transmission is close to zero, and therefore, the input for 
the GOG simulations of the Gaia RP spectra in the wavelength 
region where the transmission is high comes only from a syn- 
thetic model (the best match) that was actually used during 
the training phase of the algorithms. As a consequence, the errors 
in the parameter estimates for these stars are overly optimistic. 
The IRTF library is an intermediate case with most spectra cov- 
ering the wavelength region above 0.8 yum where most of the 
source flux is concentrated. Thus, the completion will not have 
a relevant impact on the resulting simulated RP spectrum. 

Figures 9 and 10 show the predictions of the kNN and GP 
models for the empirical libraries described in previous para- 
graphs. The predicted values of T^ff are compared with the effec- 
tive temperature assigned by the SLC cahbration to the spectral 
type provided by the empirical libraries. The model used in all 
panels corresponds to the one trained with noisy spectra of G=20 
and 70 transits. The upper panels describe the performance of 
the model when applied to noisy spectra corresponding to G=15 
and 70 transits whilst the lower panel shows the results for the 
same model appUed to GOG simulations at G=20 and 70 tran- 
sits. For each spectrum simulated by GOG, we constructed ten 
noisy replicates using the Gaia error model currently available, 
and predict the values of Tf.g and \og{g) for each one of them. 
The scatter in the predictions for these noisy replicates is better 
visible in the lower panels (spectra simulated at G=20). 

Table 3 collects the root mean squared errors (RMSE) of 
the Jeff estimates (expressed in Kelvin) obtained by the GP and 



kNN models trained with a collection of synthetic spectra of the 
BT-Settl hbrary simulated for an apparent magnitude G=20. The 
models are then apphed to the four libraries of empirical spec- 
tra simulated at 28 and 70 transits, and for apparent magnitudes 
G=15, 18, and 20. The prefix PCA refers to the models built for 
the input space of principal componets. The values in parenthesis 
correspond to the RMSE after correcting for the systematic bias 
in GP predictions. We also include the RMSE of the Bayesian 
inference modules built from the BT-Settl and COND/DUSTY 
model libraries for comparison. 

4.1. Bayesian inference 

Figure 11 shows the T^g predictions obtained from the elHp- 
soidal samplings, using the neural network trained with the BT- 
Settl grid of models and the current error model for the Gaia 
RP spectra. The predictions are obtained as before, for ten repli- 
cates of each GOG simulation of an empirical spectrum. It shows 
a tendency to predict effective temperatures around 1700 K for 
stars with spectral type L. Furthermore, there are indications of 
a second attractor slighdy above 1000 K for stars and brown 
dwarfs with spectral types between L5 and T5. These systematic 
biases can also be recognised in Figs. 8 (representing the T^fi es- 
timates from;^f^ fits) and 9 (representing the T^s estimates from 
A:-nearest neighbours), but not in Fig. 10 which corresponds to 
the GP model. The fact that the Bayesian inference shows the 
same kind of biases as nearest neighbours and the x^ fits is not 
surprising because Bayesian inference is equivalent to maximum 
likelihood estimation under flat priors such as those used in this 
work. Minimum fits and nearest neighbours are special cases 
of maximum likelihood estimation. 

To understand the nature of this trend to concentrate the pre- 
dictions by the Bayes module around 1700-1800 K, we plot in 
Fig. 12 the best fits produced by the Bayesian inference mod- 
ule (blue), the model corresponding to the spectral type assigned 
in the spectral compilations (orange), and finally, one of the ten 
noisy replicates of the spectrum simulated with GOG for G=20 
(black), for two stars with spectral types LI (left) and TO (right), 
both of which are predicted to have temperatures between 1600 
and 1700 K. For these two stars we also plot the original spectra 
for reference (bottom row). 

The reason for this systematic effect seems related to the fact 
that between 1600 and 1800 K, the model spectra undergo, es- 
pecially at log(g)=5.0, rapid changes with temperature, whilst 
before and after this range, we find plateaux where the spec- 
tra have only a mild dependence with temperature (see Fig. 13). 
From the sampUng perspective, relatively strong changes in the 
proposed temperatures falling within the plateaux regions result 
in small changes in the likelihood. Sampling in the 1600-1800 K 
region on the contrary provides a more varied range of models 
susceptible to better fit the observed noisy spectrum. 

The detailed analysis of the predictions reveals a large vari- 
ance within some of the blocks of ten noisy replicates of a 
given spectrum. A good example of this kind of problems is 
provided by SDSS J010752.33+004156.1, where six estimates 
cluster around 1036 K whilst the remaining four cluster around 
1664 K. As a reference, the spectral type quoted for this object 
is L5.5 or, equivalently, 1554 K according to the SLC calibra- 
tion. This apparent inconsistency is due to the multimodal poste- 
rior distribution, which is closely related to the likelihood land- 
scape for flat priors and noisy spectra. Figure 14(a) shows the 
log-likelihood distribution (derived from one noisy GOG simu- 
lation of the spectrum in the Leggett compilation) as a function 
of file parameters for SDSS J010752.33+004156.1. For compar- 
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Fig. 8: Effective temperatures derived from the best fits to BT-Settl models as a function of the spectral type assigned in the 
literature. Black circles correspond to the compilation by Leggett, red circles to the compilation of Keck LRIS spectra by Reid, 
orange circles to the NIRSPEC compilation, and blue ones to the IRTF compilation. The grey continuous line shows the Teff- 
spectral type calibration by Stephens et al. (2009) from optical and infrared spectra. The dashed lines represent the same calibration 
+ 250 K. The right-hand panel shows the residuals (reff(predicted)-7^eff(SLC)) with respect to the calibration. 



Table 2: Average bias (p) and standard deviation of the;^'^ effective temperature (K) fits. 



Library 


fi (BT-Settl) 


o- (BT-Settl) 


RMSE (BT-Settl) 


fi (AMES) 


o- (AMES) 


RMSE (AMES) 


Leggett 


42 


196 


199 


156 


197 


250 


Reid 


41 


72 


143 


-84 


116 


235 


NIRSPEC 


67 


177 


256 


63 


175 


251 


IRTF 


53 


63 


126 


61 


70 


140 



Table 3: Root mean square errors (RMSE) of the GP and kNN models for the prediction of Teff, appHed to the compilations of 
empirical spectra. 



Number of transits 


28 


70 


G 


15 


18 


20 


15 


18 


20 


GP-G20 


257 (192) 


260 (194) 


266 (201) 


256 (191) 


257 (192) 


260 (196) 


PCA-GP-G20 


266 (199) 


270 (203) 


281 (216) 


264 (197) 


266 (200) 


272 (207) 


kNN-G20 


209 


210 


213 


207 


209 


213 


PCA-kNN-G20 


211 


213 


215 


209 


207 


210 


Bayes BT-Settl 


230.5 


235.7 


239.0 


243.4 


241.6 


239.7 


Bayes COND/DUSTY 


252.6 


252.3 


255.0 


257.5 


257.6 


258.0 



ison, we show the equivalent plot for the model in the BT-Settl 
grid corresponding to 6> = (1550/:, 5.0) (Fig. 14(b)). Whilst the 
log-likelihood landscape for the BT-Settl model is unimodal, 
the log-likelihood landscape for the noisy spectiTim of SDSS 
J010752. 334-004156.1 shows maxima of comparable height at 
1036 and 1664 K (in addition to other local maxima in the range 
1000 K < Teff < 1500 K, 3.5 < log(^) < 4.0). Depending on the 
particular realisation of the noise, the ellipsoidal sampling will 
converge to different maxima. In practice, we find cases of clear 
bimodality in spectra of stars with spectral types between LI and 
T2 (reft between 1200 and 2000 K). 

This problem is easily solved by using priors based on the 
additional information provided by the Gaia astrometric mea- 
surements, as suggested in Sect. 3. For this case. Fig. 15 shows 
a particular choice of the prior that would discard all local max- 
ima that are inconsistent with the Gaia photometry/astrometry, 
under the assumption of negligible circumstellar and interstel- 
lar extinction. In it, we represent the distribution of absolute 



G magnitudes of the BT-Settl models as a function of the de- 
pendent parameters T^ff and log(g). The continuous line shows 
the 1-S countour of the two-dimensional Gaussian prior (with E 
being the covariance matrix). The prior is fully defined by the 
mean fiprmr and covariance Uprmr of the Gaussian prior distri- 
bution. Given a set of astrometric, photometric and spectropho- 
tometric observations such as those simulated for G=20.0 with 
GOG for SDSS J010752.33+004156.1, and the photometric and 
astrometric errors shown in Fig. 5, we calculate an absolute G 
magnitude equal to 19.8492;^pp}g4, where the uncertainties are 
derived from the values shown in Fig. 5 multiplied by five. We 
used five times the nominal uncertainties in the apparent G mag- 
nitude and 7T (the parallax) to account for the potential mismatch 
of the model-predicted magnitudes with respect to the real dis- 
tribution. We did not take the Lutz-Kelker bias into account in 
the computation of the uncertainties of the absolute G magni- 
tude. Models within the uncertainties in the absolute G magni- 
tude (19.8492+°"[^^) are shown in Fig. 15 as black circles. We 
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Fig. 9: reft predictions obtained by the kNN algorithm based on the BT-Settl grid of models for the four libraries of ground-based 
spectra listed in section 4 for G=15 (top row) and G=20 (bottom row). The colour code is the same as used in Fig. 8. The jc-axis 
shows the spectral types gathered from the literature cited in the spectral compilations. The grey line shows the Teft-spectral type 
calibration by Stephens et al. (2009) from optical and infrarred spectra (±250K, grey dashed lines). The right-hand panels show the 
residuals with respect to the caUbration as in Fig. 8. 



used the Tgft and \og(g) parameters of these models to propose 
the values of fiprior and I.prior that were used to draw the 1-E 
iso-countour in Fig. 15. This prior is eleven orders of magnitude 
larger at the 1664 K maximum than at the 1036 K one, and thus 
renders this local maximum insignificant in the posterior proba- 
bility density distribution. There is a naiTow local maximum at 
TeffSilSOO K and log(g)a;4.9 in Fig. 14 where the prior is only 
one order of magnitude smaller than that corresponding to the 
1664 K, but this maximum is never significantly sampled in any 
of our ten replicates due to its narrowness. 



4.2. Model selection 

Ellipsoidal sampling allows estimating the evidence given an ob- 
served spectrum and a model choice (COND-i-DUSTY or BT- 
Settl), as defined in Eq. 2. Therefore, it is possible to compare 
these two model families from a Bayesian perspective. Mq de- 
notes the BT-Settl model librai-y and Mi the COND+DUSTY 
library. If these were the only two alternatives, then the Bayes 
factor (BF) in support of the Mq model would be defined as the 
ratio of the respective marginal densities (evidences) of the data 
for the two models. 



BF : 



mo(s) 



(6) 



mi(s) 

If ttq and tti denote the respective prior probabilities (in our 
case TTo - TTi - I /2), the posterior probability of Mq is given by 



Pr(Mo|s) 



jtqBF 

TTqBF + TTl 



We obtain that 73.4 % of the spectra support the BT-Settl 
model librar-y against 26.6 % supporting the COND-DUSTY 
combination. For this reason we have exemplified the regression 
results with figures and discussion related to the models obtained 
with the BT-Settl library. 

4.3. Estimates of the gravity 

Figure 16 shows the T^ff-logig) diagrams obtained from the 
three models discussed in this section. The ellipses represent the 
covariance estimated from the ten noisy replicates at G=15 of 
each GOG simulations of the empirical spectra. Since we do 
not have a compilation of surface gravities available to assess 
the overall validity of the predictions, we compare the regres- 
sion values with those obtained from the fitting to the full 
resolution spectra. Figure 17 shows this comparison. The fit 
values are jittered with a Gaussian distribution of standard devi- 
ation equal to 0.2 to enhance the visibility. If the^^ estimates are 
taken as targets, only the GP model can be used to obtain very 
rough estimates of the gravity and to tag low-gravity candidates. 

The only star in our samples with an indication of low 
gravity in the comments section of the Dwarf Archives^ is the 
NIRSPEC target star 2MASS J1726000+153819. For this star, 
both the kNN model and the Bayesian estimate agree to assign a 
value of log(g)=3.5, while the GP model assigns a higher value 
log(g)=4.1. 



(7) 



* http : //spider . ipac . caltech. edu/staff/davy/ ARCHIVE/ 
index. shtml 
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Fig. 10: Teff predictions obtained by the GP model based on the BT-Settl grid of models for the four libraries of ground-based spectra 
listed in section 4 for G=15 (top row) and G=20 (bottom row). The colour code is the same as used in Fig. 8. The x-axis shows the 
spectral types gathered from the literature cited in the spectral compilations. The grey line shows the Teff-spectral type calibration 
by Stephens et al. (2009) from optical and infrared spectra (±250K, grey dashed lines). The right-hand panels show the residuals 
with respect to the caUbration as in Fig. 8. 



5. Conclusions 

We have presented the module that will be in charge of detect- 
ing and characterising ultra-cool dwarfs in the Gaia database. 
The module is subject to change and improvement, but this im- 
plementation provides the baseline performance that can be ex- 
pected from it. 

We used the current instrument models and the estimated 
spatial densities by Caballero et al. (2008) to predict the expected 
number of ultra-cool dwarfs per spectral type bin. We found that 
Gaia will be able to detect significant numbers (around or above 
ten detections) for UCDs of spectral types hotter than L6-7 v. We 
also used the BT-Settl library of synthetic spectra to define se- 
lection criteria for the UCD module such that no UCD is missed 
due to measurement eiTors. Given the spatial densities estimated 
by Caballero et al. (2008), we derived contamination rates from 
stars hotter than 2500 K in the resulting samples. 

We conducted an extensive study to find the best statisti- 
cal regression model of the relationship between the observed 
Gaia RP spectrum and the source physical parameters (T^ff and 
log(g)). We evaluated several alternatives in view of their inter- 
nal and external errors. The internal errors were estimated with 
cross validation experiments with a dataset interpolated from the 
nominal grid of models provided by the libraries of synthetic 
spectra. The external validation was carried out by applying the 
regression models to an independent set of UCD spectra ob- 
served from the ground. All these experiments were carried out 
on GOG simulations of the full-resolution spectra for a number 
of apparent G magnitudes and numbers of transits. 



As a result, we found that the expected end-of-mission error 
of the UCD module for the faintest detectable UCDs (G=20) is 
210 K for the kNN module and 260 for the GP module (207 K 
if a bias correction is applied). These performances are approx- 
imately constant as a function of G, at least down to G=15 (a 
typical apparent magnitude for the brightest UCDs in the Gaia 
catalogue), and are remarkably close to the performance of a 
simple fit with the full-resolution spectrum. 

The Bayesian inference of the source parameters shows sys- 
tematic deviations in the distribution of predicted temperatures, 
which are also apparent (although less conspicuous) in the 
fits and kNN predictions. It is also severely affected by the mul- 
timodality of the likelihood maps. The application of physical 
priors and advanced sampling techniques capable of identifying 
multiple modes in the posterior will be the subject of a forth- 
coming paper in this series. 

The logg predictions are characterised by a typical RMSE of 
0.2 dex for the GP module (0.7 dex for the kNN module) as mea- 
sured by the cross-validation experiments. Unfortunately, these 
error estimates are overly optimistic because they are derived 
from testing the regression modules on synthetic spectra and not 
on observed spectra of real UCDs. The \og(g) estimates for the 
empirical spectra observed from the ground, although broadly 
consistent with the expected distribution of values in the sam- 
ples of empirical spectra, prove that the errors quoted above are 
indeed extremely optimistic. The GP module will undoubtedly 
benefit from a more realistic (i.e. wavelength-dependent) treat- 
ment of the noise parameter. 
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Fig. 1 1 : Teff predictions obtained from the ellipsoidal samplings (BT-Settl) for the four empirical libraries of ground-based spectra. 
The colour code is the same as used in Fig. 8. The .v-axis shows the spectral types gathered from the literature cited in the spectral 
compilations. The grey line shows the Teff-spectral type calibration by Stephens et al. (2009) from optical and infrared spectra 
(+250K, grey dashed lines). The right-hand panels show the residuals with respect to the calibration as in Fig. 8. The top row 
corresponds to the spectra simulated at G=15 and the bottom row corresponds to the G=20 replicates. 
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Fig. 12: (Top panel, black continuous lines) GOG simulations of spectra contained in the Legget compilation corresponding to 
2MASS0345H-25 (left) and SDSS1511+06 (right) for G=20 and 70 transits. The blue line represents in both panels the model 
corresponding to the mode of the posterior probability as derived using ellipsoidal sampling. In orange, the model that corresponds 
to the effective temperature derived from the spectral type (LI and TO respectively) using the SLC calibration. The lower row shows 
the original spectra completed with the best model. 
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Fig. 13: BT-Settl model spectra between 1200 and 2000 K (the range of temperatures where the Bayesian module shows a clear 
tendency to concentrate predictions around 1600-1800 K) and log(g)=5.0. Each panel shows spectra in an intei-val of temperatures 
of 200 K. The colour code represents the increment in effective temperature with respect to the lowest temperature covered in the 
panel. Blue lines coiTespond to this lowest temperature and red continuous lines correspond to the lowest temperature plus 200 K. 
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Fig. 14: Log-likelihood landscapes for a G=20 noisy replicate of the spectrum of SDSS0107 (left), and the BT-Settl model for 1550 
K and log(^)=5.0 (right). 
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Fig. 15: Absolute G magnitudes for the BT-Settl model library as a function of the physical parameters T^ff and log(g) (see colour 
code at the right-hand side of the scatter plot). Models within Gabs = 19.8492+[J{J[^^ (corresponding to SDSS J010752. 33-1-004156.1 
simulated at G=20.0) are shown as black circles. The ellipse shows the 1-Z iso-contour of the Gaussian physical prior described in 
the text. 




Fig. 16: Teff-logig) predictions for the empirical spectral libraries obtained with the GP model (left), the kNN model (middle), and 
the Bayesian inference (light). The colour code is the same as in Fig. 8. The ellipses correspond to the covariance estimated from 
ten noisy replicates of the GOG simulated spectra (G=15). 
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Fig. 17: log(^) predictions for the empirical spectral libraries obtained with the GP model (left), the kNN model (middle), and the 
Bayesian inference (right), as a function of the log(§) value of the minimum fit (jittered with a Gaussian distribution of cr = 0.2). 
The colour code is the same as in Fig. 8. 
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